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Abstract 

These are notes prepared for presentation at the workshop "Chahenges in 
Granular Matter" at the Abdus Salam Institute for Theoretical Physics, Tri- 
este, August 2001. Revisions and figures will be added at a later date. 

Many features of real granular fiuids under rapid flow are replicated by a 
system of smooth hard spheres with inelastic collisions. For such a system, 
it is tempting to apply standard methods of kinetic theory and hydrody- 
namics to calculate properties of interest. The domain of validity for such 
methods is a priori uncertain due to the inelasticity, but recent systemmatic 
studies continue to support the utility of kinetic theory and hydrodynamics 
as both qualitative and quantitative descriptions for many physical states. 
The current status of kinetic theory and hydrodynamics is reviewed and in- 
terpreted, with the optimistic conclusion that much existing phenomenology 
can be placed on firm grounds for further understanding of complex states 
observed in rapid granular fiow. 

I. INTRODUCTION 

Granular media in rapid, dilute flow exhibit a surprising similarity to ordinary fiuids and 
the utility of a hydrodynamic description for such conditions has been recognized for many 
years [||]. This phenomenology has come under scrutiny in recent years with questions about 
the domain of its validity and the associated constitutive equations appearing in the hydro- 
dynamic equations Answers to such questions can be found in a more fundamental 
microscopic description where the tools of nonequilibrium statistical mechanics are available 
for critical analysis. An intermediate mesoscopic description between statistical mechanics 
and hydrodynamics is that of kinetic theory, whose applicability to granular matter also 
poses questions. Attention here is restricted to an idealized microscopic model system that 
captures only the primary feature of grains: inelastic collisions. The system is composed of 
smooth, hard spheres whose energy loss on binary collisions is characterized by a restitution 
coefficient < a < 1. This single parameter distinguishes the ideal granular fluid (a < 1) 
from the ideal normal fluid (a = 1). For normal fluids, extensive studies during the past 
thirty years confirm the validity of a wide range of many-body methods to analyze states 



near and far from equilibrium. These studies include applicability of kinetic theory (Boltz- 
mann, Enskog, and corrections due to correlated collisions), and on longer space and time 
scales, a hydrodynamic description. The foundations for kinetic theory and hydrodynam- 
ics, including practical calculation of properties such as the equation of state and transport 
coefficients, is most complete for the elastic hard sphere system. 

How much of this body of knowledge can be translated to rapid flow granular media 
using the same system with the only change being a < 1? At the most fundamental level of 
Newton's equations there is little evident change for the fluidized state; the microstate of the 
system at time t is completely characterized by the positions and velocities of all spheres. The 
sequence of free streaming and instantaneous velocity changes on binary collisions determines 
the evolution uniquely for given initial conditions at t' < t. Only the binary collision rule 
changes, such that the relative velocity gj^ = Vj — Vj for particles i and j after the collision 
is given by 



where a is a unit vector along the line from particle j to i. This is the usual rule for reflection 
of relative velocities except that the magnitude of the normal component \gij ■ a-\ after the 
collision has been decreased by a factor a to a \gij ■ a\. Thus, the standard event driven 
molecular dynamics (MD) simulation method developed for normal fluids can be applied 
directly to granular fluids with this small change in the collision ruleQ. Given the microstates 
and dynamics prescribed in this way, the tools of nonequilibrium statistical mechanics can 
be applied as well. In particular, the formulation of a kinetic theory and its application to 
derive the equations of fluid dynamics for granular systems can be addressed. The objective 
here is to summarize the current status of such investigations and to address some of the 
early and continuing concerns about the validity of kinetic theory and hydrodynamics for 
granular media. Relevant literature prior to 1990 has been reviewed by Campbell ||^. More 
recent reviews addressing the statistical mechanical foundations have been given by Ernst 
and by van Noije and Ernst 0. The reader is advised to consult the extensive references 
in these articles. Here, a personal perspective on the current status of kinetic theory and 
hydrodynamics for granular gases is provided as a supplemental orientation to this rapidly 
growing fleld. The references given are selective and apologies are offered at the outset to 
the many other important contributions not recognized explicitly. 

Contributions to this fleld come from increasingly diverse areas: chemical engineering, 
fluid dynamics, statistical mechanics, condensed matter physics, computational physics. 
Each has its priority for problems to be solved, in some cases focused more on urgent prac- 
tical concerns while in others addressing more abstract "points of principle". The material 
selected for discussion here reflects the author's personal view that a fundamental (micro- 
scopic) basis for rapid flow granular media can be found within the general framework of 



The phenomenon of "inelastic collapse" can stop the MD simulation due to condensation of small 
arrays of particles with rapid collisions at high density and strong dissipation ||2| . This is an artifact 
of the instantaneous hard sphere collision that appears to have no macroscopic consequences, and 
which can be avoided in several different ways. It will not be considered further in the discussion 
here. 
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nonequilibrium statistical mechanics, and that its exphcation will impact practical applica- 
tions under realistic conditions. This perspective approaches granular media via idealized 
models for which controlled molecular dynamics simulation, many-body theory, and real ex- 
periments reasonably can be performed to answer specific questions with clear consequences. 



II. KINETIC THEORY AT LOW DENSITY 

It is taken for granted in all of the following that a microscopic formulation of non- 
equilibrium statistical mechanics for inelastic, smooth hard spheres is justified Con- 
cerns at this level are the same as those for open systems (inherently non-equilibrium states) 
with elastic collisions - e.g. the required "mixing" property for the dynamics in phase space 
to justify the relevance of average values and fluctuations about them. Once it is recog- 
nized that granular states are always "non-equilibrium" the distinction between systems 
with elastic and inelastic collisions diminishes considerably, as emphasized below. A first 
order question is whether the reduction of non-equilibrium statistical mechanics to a more 
practical kinetic theory description can be justified. This is a difficult problem even for the 
case of elastic collisions which, as noted above, has been addressed extensively and success- 
fully for the hard sphere system over the past thirty years. It is not the purpose here to 
summarize the status of that research but rather to explore any additional questions and 
concerns for its application to granular gases arising from the inelasticity of the collisions. 
The general conclusion below is that while there are significant differences in many derived 
properties of interest, there are few new qualitative limitations on the equations describing 
these properties based on analytic and simulation studies to date. Some of the confusion and 
disagreement in the literature is due to differing perceptions of the necessary conditions for 
a kinetic theory (e.g., Boltzmann or Enskog kinetic equations), and the necessary conditions 
for a hydrodynamic description (e.g., Navier-Stokes or non-linear rheological equations). In 
principle and in practice, the validity of kinetic theory and the validity of hydrodynamics 
are independent questions. However, for gases (no compact structures) it is expected that 
existence of a kinetic theory should provide the basis for investigating hydrodynamics as 
well. 

To set the stage for the discussion here, a brief formal "derivation" of the kinetic de- 
scription for both the distribution function and fluctuations at low density is given as an 
illustrative example. A small parameter expansion of the reduced distribution function leads 
to a formal solution to the BBGKY hierarchy P,pi| . The analysis is similar to an expansion 



proposed by Grad for the hard sphere gas with elastic collisions |Tl|]. In this approach, there 
is no reference to concepts such as "approach to equilibrium", Maxwellian distribution, or 
"molecular chaos", and the distinction between inelastic or elastic collisions plays no explicit 
role. Thus, superflcially at least, it appears the basis for the low density kinetic theory is the 
same in both cases. Closely related results follow from a cluster expansion of the reduced 
distribution functions, again with no formal reference to whether the collisions are elastic or 
inelastic i,|,|12[. 

The s-particle reduced distribution functions, f^^\xi, ■ ■ ■ ,Xs,t), obey the BBGKY hier- 
archy where Xi = (q^, Vj). A dimensionless form of this hierarchy is obtained by scaling the 
space and time with the mean free path i = l/{na'^) and the mean free time to = ^/vq. 
Here, n is the density, cr is the hard sphere diameter, and Vq is some characteristic veloc- 
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ity. Similarly, the reduced distribution functions are scaled with (n/fg)*. The resulting 
dimensionless BBGKY hierarchy has the form 

^ s ' _ \ 

a* + E ■ V. - ^ Til, j) / W (xi, • ■ ■ , a;,, t) 

i=l i<j J 

= ^2 [ dxs+iT{i,s + l)f^'^'\xu - ■ ■ ,Xs+i,t). (2.1) 
i=i 

The binary scattering operator T{i,j) for a pair of particles {i,j} is defined by 

T{t,j) = Jda 0(g,,-?)(&,-5)[a-25(qy - \a)br^' - 5(q,, + \a)], (2.2) 

where a is a unit vector along = q^ — q^ and da denotes a two dimensional solid angle 
integration over the sphere for particles at contact. Also, gij = Vj — Vj is the relative velocity, 
and b~j^ is the scattering operator defined for any function X(vj, v^) by 

6,r.ix(v„v,)=X(v:,v;.). (2.3) 



The "restituting" velocities (i.e., the inverse of those in (|1 . 1| ) ) are 



Vi = Vi-^—{g,j-a)a, Y.=Vj + ^—{g,j-a)a. (2.4) 

The a dependence of T{i,j) contains all aspects of the inelasticity, and plays no explicit role 
in the following expansion. 

In this dimensionless form the BBGKY hierarchy depends on the single dimensionless 
parameter X = a/i = na^, the ratio of the "force range" to the mean free path. This 
parameter is small at low density, suggesting an expansion for a solution to the entire 
hierarchy as a power series in A. The dependence on A occurs explicitly as shown on the left 
side of ( |2.1| ) and implicitly through the finite separation of the colliding particles in T{i,j). 
The structural features of the expansion in A are simplest if it is performed at fixed T{i,j). 
The s-particle reduced distribution functions are taken to have the representation 

/(^Hxi, ■ ■ ■ , x„ t) = fi'\x,, ■■■,Xs,t) + \^f['\x^, ■ ■ ■ , x„ t) + ... (2.5) 

It is then readily shown that the hierarchy is solved exactly to order A^ in the form 

#(a;i,---,a;„t) = n/o^')(a;,,t), (2.6) 

i=l 



/r^(xi,---,x„t)=^n/ona:.,t)/r(x„t)+5: n fr{xk,t)G{x,,x,,t), (2.7) 

where the expression for fi^'' holds for s > 2. Thus, the reduced distribution functions for 
any number of particles is determined as a sum of products of the single particle functions 
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fQ^\xi,t) and fi^\xi,t), and the pair function G{xi,X2,t). These are determined from the 
set of three kinetic equations 

+ vi ■ Vi] tf)(xi,t) = J{x„t\ /o«). (2.8) 



(^^ + vi • Vi - Ji + V2 ■ V2 - G{x,,X2, t) = T(l, 2)f^'\xu t)g\x2, t) (2.9) 

(^^ + vi ■ Vi - fi'\x,,t) = J dy2T{l,2)G{xi,X2,t), (2.10) 

Here J{xi,t \ /o^'') is the Boltzmann-Bogohubov colhsion operator and Ji, defined over 
functions of Xi, is its hnearized form 

J{xi,t\ tf^) = / dx2T{l,2)f^'\x^,t)f!,'\x2,t) (2.11) 



hh{x^) = f dx2T{l,2) {jf\x,,t)h{x2) + h{x^)ff\x2,t)) . (2.12) 



These low density results ( |2.8| )-( PTTI1| ) are remarkably rich. The leading order distribution 
function Jq^^ is the solution to the Boltzmann equation. The two particle correlations are 
generated from the uncorrelated product of Boltzmann solutions through inelastic binary 
collisions T{l,2)fQ^\xi,t)fQ^\x2,t) on the right side of (P^). Finally, corrections to the 
Boltzmann solution due to correlations are given by a coupling of the distribution function 
to the correlations in (|2.10| ) (the so-called "ring" recollision effects). In the following sections, 



it is noted that the solutions and implications of these kinetic equations can be quite different 
for elastic and inelastic collisions. But these differences come from the equations themselves 
and should not be interpreted as signatures of their failure to apply. For example, at a = 1 a 
possible solution for an isolated system is fQ^\xi,t) — > fM{vi), G{xi,X2,t) = = f[^\xi,t), 
where /m is the Maxwell-Boltzmann distribution. Equation ( ^.1(J| ) supports because 
energy is conserved, and G{xi,X2,t) = because T(l, 2)/j\/(vi)/j\/(v2) = for the same 
reason. Since energy conservation no longer holds with a < 1 it is not surprising that the 
isolated system does not approach equilibrium, the Maxwellian is not a stationary solution, 
and that finite correlations exist. Indeed, the extent to which such predicted differences 
agree with observations from molecular dynamics provide support for the kinetic theory, not 
limitations on it as is sometimes implied. 

Clearly, the above derivation has not restricted this kinetic description to isolated sys- 
tems. In fact, the most interesting cases of practical interest are response to boundary 
conditions and/or external fields. The similarities between normal and granular fluids is 
closest for such " nonequilibrium" conditions. Too often, properties of granular gases are 
contrasted only to those of the equilibrium state for normal gases. It is important to note 
that practical access to the solutions to the above kinetic equations is possible for a wide 
range of conditions by direct simulation Monte Carlo (DSMC) ||13|. Analytic results for 



5 



complex states have been obtained using accurate kinetic models for the above equations 
as discussed in Section 7. The benchmark for validity of these kinetic equations is 



molecular dynamics (MD) simulation of the idealized granular fluid microdynamics. Fre- 
quently, these simulations are performed at higher densities than can be supported by the 
low density approximation. However, there is a corresponding Enskog approximation that 
extends the above low density kinetic theory to moderate densities as well [[7|,|15|. It is appro- 
priate to emphasize also that the terminology "kinetic theory" is not limited to Boltzmann 
and Enskog equations. Indeed, at high densities and long times it is known that there are 
important (even dominant) corrections due to correlated many-body collisions of the type 
described by ( |2.10| ) at low density. Thus, observed deviations from Boltzmann-Enskog for 
a < 1 do not necessarily imply a failure of a kinetic theory description in general. 



III. HOMOGENEOUS COOLING STATE 

An isolated homogenous fluid with elastic collisions rapidly approaches the stationary 
equilibrium Gibbs state after a few collisions. The single particle distribution function is 
Maxwellian and there are no correlations on length scales large compared to the hard sphere 
diameter. In contrast, for inelastic collisions no stationary state is possible for the isolated 
system. Instead, it is straightforward to show that the temperature (defined in the usual way 
from the average kinetic energy) is monotonically decreasing due to the inelastic collisions. 
Still, on the same time scale as the approach to equilibrium for elastic collisions, the granular 
fluid approaches a special "homogeneous cooling state" (HCS) for which all time dependence 
of the distributions occurs only through this cooling temperature. This implies a scaling 
form for the single particle distribution 

f!,'\v,t) = v^'it)n<Piv/voit)), vl{t)=2T{t)/m. (3.1) 

The pair correlations have a similar scaling form. The existence of a scaling solution and 
corresponding time dependence of T(t), its approach from general homogeneous states, the 
form of (f ) , and the space dependence of correlations in this state provide specific tests 
for the kinetic theory. The scaling property implies an exact time dependence T{t) = 
T(0)(1 + ct)~^ (Haff's law) which has been confirmed by MD ]16[ and by DSMC solution 
to the Boltzmann and Enskog equations [0. In both the MD and DSMC the approach to 



a scaling solution is observed to occur rapidly after a few collisions, similar to the approach 
to equilibrium for normal fluids. The velocity dependence of (p ('^) has been calculated 
approximately from the the Boltzmann and Enskog equations by an expansion about the 
Gaussian (a = 1 form) in polynomials p!8|JT9|| . Deviations from the Gaussian have been 
measured by MD and DSMC and compared with a low order truncation of the expansion. 
The agreement is excellent for all 0.3 < a < 1. Finally, the form of this scaling distribution 
has been analyzed from the Boltzmann equation for large velocities and predicted to have 
an exponential rather than Gaussian form |]T9|. This also has been confirmed in detail by 
DSMC and the effect has been observed by MD simulation as well . 



Recently, the scaling solutions to the Boltzmann equations for a mixture have been 
studied [^. A new effect is the occurrence of different temperatures for each species. 



although each with the same cooling rate given by Haff's law. This difference is a reflection 
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of the failure of energy equipartition in the HCS, an exact property of the Gibbs state. The 
polynomial solutions for the mixture distribution functions can be constructed in the same 
way as for the one component fluid. The rapid approach to this HCS, the dependence of the 
different temperatures on species parameters (size, concentration, mass), and the velocity 
dependence of the scaling solutions have been confirmed in detail by DSMC ||22|, for a > 0.5. 
To date there have been no MD simulations for comparison. 

Correlations in the HCS are described by Eq.(|2.9|) for G{xi,Xj,t) = /(^^(xi, X2, t) — 
f^^\xi,t)f^^\x2,t). The origin of these correlations is the source term in this equation 
T(l, 2)0 (t>i/t>o(t)) 0(f2/fo(t)) 7^ 0. This vanishes for elastic collisions because in that case 
becomes the Maxwellian, and T(12) regains the detailed balance property. For a < 1 this 
finite source gives non-zero short ranged spatial correlations ~ S{qij — cr). The generators 
for the two particle dynamics on the left side of ( p.9| ) include in their spectrum the inelastic 
fluid hydrodynamic modes (discussed below). These modes give rise to long ranged spatial 
correlations in G{xi,Xj,t), including algebraic decay at large distances for velocity moments 
corresponding to density, energy, and momentum correlations. A detailed study of the 
predictions for these correlations based on ( |2.9| ) has been given by Brey et al. |T^. These 
calculations require analysis of the spectrum for vi ■ Vi — /i at long wavelengths, and confirm 
the existence of hydrodynamic modes as spontaneous fluctuations in the granular gas. These 
modes agree with those obtained for the decay of imposed spatial gradients described in 
the Section 6, and this spectral analysis can be viewed as an independent confirmation of 
hydrodynamics (in effect, this is Onsager's regression hypothesis extended to granular fluids). 
The dispersion relations for these modes as a function of a and wavevector is more complex 
than for normal fluids and the spatial dependence cannot be obtained from simple scaling 
without specification of a domain for a. The calculated space dependence of the correlation 
functions is in excellent agreement with DSMC simulation of the kinetic equation ( |2.9| ). 

Correlations on such macroscopic scales in normal fluids also can be described by the 
usual hydrodynamic equations extended to include thermal fluctuations characterized by 
white, Gaussian noise. Averaging products of solutions to these stochastic equations over 
the noise provides the correlation functions, van Noije et al. |2^ have extended this approach 



to granular fluids as well. The agreement with the kinetic theory and DSMC results above 
is good, except for large inelasticity where the use of elastic fluid transport coefficients in 
23| is not justified. Furthermore, these results are in good agreement with MD simulations. 



In summary, the consistency of kinetic theory and hydrodynamics to predict a new effect 
for granular gases - long range correlations in an isolated homogeneous system - appears to 
be established. 

Spatial correlations on short length scales have been described in detail by Lutsko pi,p3[. 



The analysis is based on an exact boundary condition for hard spheres at contact p6| 



6 (qi2 - 0-) (o- • gi2) (cr ■ gi2) f^^^\xi, X2, t)d\idv2 
= 6 (qi2 - cr) {-a ■ g[,) 6 {-a ■ g[,) {x\, x'^, t)dv[dv', (3.2) 

The left side describes the flux of particles on the post-collision hemisphere. The right 
side is the distribution of particles on the pre-collision hemisphere as a function of their 
restituting velocities, and the equality represents their equality under a two particle collision. 
Then using the properties —a ■ = a~^cr ■ gi2 and dw'^d'w^ = a~^d\id\'2 the half space 
distributions at qi2 = a are seen to be related by 
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{(T ■ gi2) f^''\xi, X2, t) = a-%-'e {-a ■ gi2) f^'\xi, X2, t) (3.3) 

Lutsko observes that the class of generahzed mean spherical approximations (GMSA) for 
the equilibrium radial distribution g (qi — q2, a = 1) function requires primarily the value of 
this function at contact and the equation of state. Consequently the GMSA can be extended 
to the HCS for g (qi — q2, a < 1) with appropriate changes in these quantities. The value 
at contact is obtained by using the approximation 

e (-qi2 ■ gi2) f^'Kx,, X2, t) ^ e (-qi2 ■ gi2) 9 (qi - q2, 1) f^'\x,,t)f^'\x2, t) (3.4) 

on the pre-collision hemisphere, i.e. neglect of velocity correlations on this hemisphere. 
The boundary condition (|3.2|) then gives the remaining half of the distribution on the post- 
collision hemisphere with velocity correlations determined from the binary collision. The 
equation of state is determined from the Enskog kinetic equation which can be obtained in 
this context as follows . The exact first equation of the BBGKY hierarchy (Eq. ( p.l|) with 



s = 1), is written in the equivalent form 

+ vi ■ Vij f^^\xi,t) =cr^ J J da 5 [qu - (tB) {B ■ gu) a^^^y^ + 1 

xe (-a ■ gi2) /2(xi,X2,t). (3.5) 

This shows that any approximation for f2{xi, X2,t) need only be imposed for the pair of 
particles at contact and on the pre-collision hemisphere. In particular, the use of ( |3.4| ) in 
( p.5|) leads to the Enskog kinetic equation. This route to the Enskog (and hence Boltzmann 
equation at low density) provides a precise meaning to the phenomenological terminology 
"molecular chaos assumption ": neglect of two particle velocity correlations for a pair of 
particles at contact on the pre-collision hemisphere. The pressure in the HCS follows from 
the Enskog equation by evaluating the trace of the average momentum flux. In this way, 
the radial distribution function at contact and pressure are found to be g {qu = a,a) = 
^2^9 (Q'i2 = cr, a = 1) and p (a) = nksT 1 + ^^^^^nna^g (gi2 = a,l) , respectively. Using 
9 {qu = c", «) and p{a = 1) in the GMSA Lutsko obtains a model for g (gi2, a) as a function 
of gi2, a, and the density. The agreement with MD simulation is very good even at extreme 
conditions of high density, na^ = 0.5 and strong dissipation a = 0.5, for the domain tested 
0" ^ Q'i2 < 3o". This is an indirect confirmation of the molecular chaos assumption ( ^.4| ), 
supporting as well the basis for the Enskog kinetic equation. In detail, however, there 
are deviations for a < 0.7 that have a system size dependence being larger for the larger 
system (in particular, the molecular chaos value at contact appears to fail). The explanation 
for this is a long wavelength hydrodynamic shear mode becoming unstable at about a = 
0.7, as discussed in the next section. When this is compensated for in the simulation the 
system size dependence no longer occurs, but significant deviations from the "molecular 
chaos" approximation remain at high density for a < 0.6. A similar limitation on the 



"molecular chaos" assumption was observed in a two dimensional MD simulation |^ where 
the distribution of impact parameters on the pre-collision hemisphere was found to become 
non-uniform (particularly for large impact parameters) at values of the restitution coefficient 
for which significant shear fluctuation is observed. 

The validity of ( |3.4| ) has been tested by MD simulation more directly , although for the 
limited domain 0.98 < a < 1, 0.05 < na^ < 0.2. These authors study a generalization of the 
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pair correlation function obtained by integrating /^^^ over the center of mass velocity and the 
magnitude of the relative velocity. This defines a pair correlation function g {qi2 = cr, 6, a) 
depending on the angle between the relative velocity and the relative coordinate. The 
MD results for g {qu = a, 9, a) as a function 9 for na^ = 0.1 and a = 0.98 show good 
agreement with the molecular chaos approximation except near 9 = tt/2 (large impact 
parameters) . To resolve small differences the pressure is measured directly from the virial and 
an order parameter for pre-coUision velocity correlations is measured. The pressure shows 
a systematic deviation from the Enskog pressure of less than one percent for the density 
range considered. The order parameter for velocity correlations has a significant system size 
dependence and an extrapolation of MD results to large systems is required. Signatures of 
finite pre-coUision velocity correlations are found in this way, indicating deviations from ( ^.41 ). 
In spite of the good agreement found for g {qi2 = cr, 9, a), the authors emphasize limitations 
of the Enskog kinetic theory based on the pressure and order parameter measurements. To 
put this in context, it should be recognized that the same small corrections to the Enskog 
equation can occur at these densities also with elastic collisions, except for the equilibrium 
state. Indeed, they can be calculated from Eq. (|2.9|) for a = 1. Once it is recognized 
that the HCS is a nonequilibrium state the role of the Enskog kinetic theory as a good 
approximation at low and moderate densities is seen to be supported for both elastic and 
inelastic collisions. However, it appears from the discussion of the last paragraph that the 
range of densities for which this is true decreases with increased dissipation. 

The strongest objections to the above description for properties of the HCS and the 
validity of the Boltzmann equation for granular gases have been offered in reference PU 



However, both the DSMC simulation and theoretical analysis presented there have subse- 
quently been criticized as invalid PT| , P^ . 

Related studies of homogeneous steady states (HSS) also have been studied as tests of 
kinetic theory [^-3^. Such steady states are produced by insertion of a random external 



force to accelerate the particles and hence compensate for collisional cooling. This has the 
advantage of avoiding the intrinsic time dependence of the HCS, but at the price of intro- 
ducing unknown new effects induced by the external force. It should be noted that the time 
dependent HCS can be mapped formally onto a stationary state for time independent study 
as well, by a suitable change of variables for the velocities and time There are consid- 
erably more MD results for the HSS, with similar conclusions regarding the deviations of 
the velocity distribution from Maxwellian, excess population at large velocities, and velocity 
correlations that are generally consistent with kinetic theory as described above for the HCS. 
Fluctuating Navier-Stokes order hydrodynamics, including noise from the external random 
force predicts well both spatial correlations and fiuctuation-renormalization of mean values 

MM- 

This section has focused on the velocity distribution and correlations for the simplest 
state of a granular gas, isolated and uniform. Significant differences are observed from the 
corresponding equilibrium state for the same gas with elastic collisions. In some of the 
literature quoted the existence of such differences from the equilibrium state interpreted 
as a failure of conditions for kinetic theory for granular gases. As the discussion above 
shows, it is quite the other way around - these differences are quite well predicted from the 
kinetic theory. Deviations from Enskog kinetic theory at higher densities are expected due 
to correlated collisions, just as for fluids with elastic collisions There is some indication 
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that the effects of correlated colhsions are enhanced at strong dissipation, and hence appear 
at lower densities, due to unstable fluctuations in the transverse flow field as discussed in 
the next section. It would be interesting to translate this qualitative statement into a more 
quantitative one by careful comparison of such effects for normal fluids in a nonequilibrium 
state with those of the HCS, since they are certainly present in both cases. 



IV. INSTABILITIES AND CLUSTERING 

It was first observed in MD simulations that the HCS is unstable to sufficiently long 



wavelength perturbations For systems large enough to support such spontaneous fluc- 
tuations the HCS becomes inhomogeneous at long times. In MD simulations the inhomo- 
geneities may grow by the formation of clusters, ultimately aggregating to a single large 
cluster ; if cluster growth is suppressed a vortex fleld may grow to the system size where 
periodic boundary conditions can induce a transition to a state of macroscopic shear. The 
mechanism responsible for the growth of inhomogeneities can be understood at the level 
of Navier-Stokes hydrodynamics where linear stability analysis shows two shear modes and 
a "heat" mode to be unstable. Quahtatively [Q, spontaneous vortex fluctuations of suffi- 
cient size grow due to this instability and provide excess local temperature and pressure; a 
particle flux is then established from high to low shear domains increasing the density and 
cooling rate; this further enhances the pressure difference to continue the particle flow. A 
more quantitative description poses a challenging test for both kinetic theory and hydrody- 
namics. 

Perhaps the simplest signature of inhomogeneities is a deviation of the global temperature 
from the form given by Haff's law. In the presence of local flow fields the average kinetic 
energy for the system decays more slowly than that implied by the homogenous scaling form 
If the local flow flelds are assumed to be dominated by vortex flow, the decay law 



can be calculated from hydrodynamics The agreement of the slower decay predicted in 



this way with MD simulation is excellent, conflrming both the instability and the relevance 
of hydrodynamics for its description. Another signature of vortex growth is provided by 
the time dependence of the transverse velocity fleld spatial correlations. As discussed above 



they can be calculated from the kinetic equation ( p.9|) as was done in references [|^,^, or 
from fluctuating hydrodynamics as in references |2^,^. The amplitude of the correlations 
in the transverse velocity fleld is seen to increase by an order of magnitude on the time scale 
for growth of the hydrodynamic shear mode. The agreement between MD, DSMC, kinetic 
theory, and hydrodynamics is very good. 

A more complete investigation of the mechanism for cluster formation has been carried 



out by Brey and collaborators based on DSMC solution to the Boltzmann equation 
and an approximate solution to the nonlinear Navier-Stokes hydrodynamic equations de- 
rived from it (Section 6). The linear stability analysis for the hydrodynamic equations 
provides a prediction for the critical wavelength for the instability as a function of a. This 
was measured by DSMC for the energy at flxed a to determine the smallest system size for 
which secular growth occurs. The agreement between hydrodynamics and simulation for this 
critical wavelength is excellent for all 0.7 < a < 0.9, i.e. including strong dissipation. Next, 
the evolution of an initial flow fleld with sinusoidal space dependence of wavelength A was 
considered. As expected from the linear stability analysis, the transverse component grows 
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in time preserving the same imposed wavelength. In contrast the initially uniform density 
and energy becomes sinusoidal at half the original wavelength, suggesting a bilinear coupling 
of these fields with the transverse flow field. Indeed, this is consistent with the mechanism 
proposed by Goldhirsch, Tan, and Zanetti described above, where the hydrodynamic 



viscous heating provides such a bilinear coupling in the temperature equation. A detailed 
comparison between a nonlinear hydrodynamic solution and DSMC for the growth of this 
second harmonic in the density and temperature shows excellent agreement This con- 



firms that the mechanism for initial cluster growth is initiated by the unstable shear mode 
and its nonlinear coupling to the density. 

There are also MD confirmations of a hydrodynamic description for the instabilities. 
Deltour and Barrat study the onset of shear and cluster instabilities for a two dimensional 



system, showing good agreement with predictions of linear hydrodynamics. A nonlinear 
hydrodynamic analysis has been compared with MD simulation at much weaker dissipation, 
a > 0.95 The density and temperature are assumed to be slaved by the linearly unstable 



transverse velocity field and quadratic nonlinearities are retained. For fixed system size there 
is a critical value of a below which the nonlinear equations support a stable inhomogeneous 
macroscopic shearing state. The critical value of a and the values of the Fourier components 
for the flow field, density, and temperature are found to be in very good agreement with the 
MD simulation results. 

In summary, for sufficiently large systems the HCS develops spatial inhomogeneities at 
times long compared to that required for growth of a local shear fluctuation. The exis- 
tence, onset, and early evolution of the instability is well-described by Navier-Stokes order 
hydrodynamics. The kinetic equation also can be used to explore the late stage evolution as 



well via DSMC. It appears |^2[ that the density field approaches an inhomogeneous steady 
form, independent of initial conditions, for fixed a. This suggests that the isolated system 
supports a more complex, but stable, state whose form is yet to be clarified. MD simulations 
do not yet confirm this, perhaps due to the smaller system sizes and time scales considered. 



V. HYDRODYNAMICS 

Consider now a spatially inhomogeneous state, created either by initial preparation or by 
boundary conditions. In the bulk, there are local balance equations for the density n{r,t), 
temperature T{r,t) (or energy density), and flow velocity U(r, t) 

Dtn + nV - V = 0, (5.1) 



DtT + ^ {Pi,d,U, + V • q) = -TC, (5.2) 

DtUi + {mn)-^djPij = 0, (5.3) 

where Df = dt + V ■ V is the material derivative, Pij{r,t) is the pressure tensor and q(r, t) 
is the heat flux. The form of these balance equations is the same as for fluids with elastic 
collisions except for the source term on the on the right side of (|5.2|) due to the dissipative 
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collisions, where ( oc {1 — a^) is identified as the cooling rate. Haff's law follow's directly 
from (|5.2| ) for spatially homogeneous states and the scaling law ( cx y/T. The fluxes Pjj, q 
and the cooling rate ( are given as explicit low degree moments of the distribution functions 
and f(^\x,,X2,t). 

These balance equations are an exact consequence of the Liouville dynamics for the sys- 
tem. Their utility is limited without further specification of Pij, q, and ( which, in general, 
have a complex space and time dependence. However, for a fluid with elastic collisions this 
dependence " simplifies" on sufficiently large space and time scales where it is given entirely 
through a functional dependence on the fields n, T, and U. The resulting functional depen- 
dencies of Pij and q on these fields are called constitutive equations and their discovery can 
be a difficult many-body problem. The above balance equations, together with the consti- 
tutive equations, become a closed set of equations for n, T, and U are called hydrodynamic 
equations. This is the most general and abstract notion of hydrodynamics, which encom- 
passes both the Navier-Stokes form for small spatial gradients and more general forms for 
nonlinear rheological transport. The primary feature of a hydrodynamic description is the 
reduction of the description from many microscopic degrees of freedom to a set of equations 
for only five local fields. 

How does this reduction come about? For the case of elastic collisions the chosen fields 
are local densities of globally conserved number, energy, and momentum. Consequently, the 
dynamics of these variables can be made as slow as desired by considering long wavelength 
phenomena. In this way a time scale can be chosen such that all other excitations have 
decayed to zero, leaving a time dependence for all properties only through these conserved 
densities. In practice, the hydrodynamic fields dominate for times large compared to the 
mean free time, and for wavelengths large compared to the mean free time. This explains 
the wide domain of applicability for a hydrodynamic description of physical phenomena. 

The extension of these ideas to granular flow raises several questions: 1) since energy 
is not conserved, why should the equation for T be included in the set of slow hydrody- 
namic fields?, 2) is the new time scale associated with the cooling rate ( microscopic or 
macroscopic?, 3) without any equilibrium Gibbs state, what is the reference state toward 
which these fields are relaxing? These issues are quite different from the those regarding the 
validity of a kinetic theory, so to address them consider values of the density and restitution 
coefficient for which the Enskog kinetic equation is reliable. This is given by the exact first 
BBGKY hierarchy equation (|3.5|) with the molecular chaos assumption 

^ + vi ■ Vi^ /(n, vi, t, ) = JeItu vi|/(t)], (5.4) 

where Je is the Enskog collision operator, 

JE[n, vi|/(t)] = a^J rfv2 J da Q{a ■g){a- g) {a-^f^'\r,, n - a, y[, V„ t) 

-/(2)(ri,ri + o-,Vi,V2,t)}, (5.5) 

/(^)(ri, ra, Vi, Vg, t) = g{ri, r2\n{t))f{ri,\i,t, )/(r2, Vg, t). (5.6) 

(Here and in the following the superscript on /^^^ will be suppressed to simplify the notation). 
The balance equations in the form j ^.lf) - (^.3|) above follow directly from this equation 
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without further approximation. Similarly, Pij, q, and ( are found to be bilinear functionals 
of / 

Pij = J d^mViVj f{r,Y,t) + ^^"^^^Z ^^vi J dva J dQ 6(5- ■ g)((5- ■ gfaidj 

X /'rfA/(2)(r-(l-A)o-,r + Acr,vi,V2,t), (5.7) 
Jo 

q = J dv^mV^Vf{r,Y,t) + ^^ma' J dv, J dw^ j dQ ■ g){{a ■ g)\G ■ 

X /'dA/(2)(r-(l-A)o-,r + A<T,vi,V2,t), (5.8) 

JO 

C=(1-«')y^/ rfvi I rfv2 |rffie(^-g)(a-g)3/(2)(r,r + ^,Vi,V2,t), (5.9) 



where V = v — U(r, t) is the peculiar velocity and G = |(Vi + V2). Equations ( ^.1|) - 
( p.9| ) provide an unambiguous basis to study the existence and criteria for a hydrodynamic 
description. Although these equations are complicated, DSMC provides direct numerical 
access to solutions to the kinetic equation to determine the fields and fluxes without any 
assumptions regarding hydrodynamics for critical tests and benchmarks. 

Hydrodynamics results from the balance equations supplemented with constitutive equa- 
tions for Pij, q, and C,. Since the former are exact, the critical issue is the existence and 
form of the constitutive equations. It is clear from (|5.7| ) - ( ^.9|) that they will be obtained if 
the Enskog equation admits a "normal" solution, whose space and time dependence occurs 
entirely through its functional dependence on the fields 

/(r,v,t) = F(v|n,T,U) (5.10) 

The fluxes and cooling rate then inherit this space and time dependence and become con- 
stitutive equations. The space and time dependence of the fields follows from solution to 
the resulting hydrodynamic equations to complete the self-consistent description of F. A 
hydrodynamic description for granular gases can be justified if it can be shown that such 
a "normal" solution exists, and that it represents a wide class of more general solutions on 
a large space and time scale. This is a difficult task even for the case of elastic collisions, 
and detailed results exist only for a few special nonequilibrium states using model kinetic 
equations (see Section 7 below). For the Boltzmann and Enskog equations partial answers 
exist for gases with elastic collisions and states with small spatial gradients. In that case. 



the Chapman-Enskog method generates the normal solution explicitly by a systematic 
expansion in the small gradients . To lowest order the resulting constitutive equation for the 
pressure tensor is Newton's viscosity law, while that for the heat flux is Fourier's law. The 
pressure and associated transport coefficients in these expressions (shear and bulk viscosi- 
ties, thermal conductivity) also are given by the method in terms of the Enskog collision 
operator. The balance equations become the familiar Navier-Stokes hydrodynamic equa- 
tions. The dominance of the hydrodynamic description at large space and time scales is 
justified by a study of the spectrum of the linearized Enskog equation. It can be shown that 
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there are five smallest eigenvalues (hydrodynamic modes) that scale with the wavelength 
such that they are isolated from all other eigenvalues |^5| for wavelengths long compared 
to the mean free path. Thus a general excitation of the system will be dominated by the 
hydrodynamic modes after a short transient time during which all other excitations decay. 
It is verified that the spectral hydrodynamic modes are the same as those obtained from 
the linearized Navier-Stokes equations obtained by the Chapman-Enskog method, at long 
wavelengths. 

While questions of convergence of the Chapman-Enskog solution and the nature of the 
constitutive equations and boundary conditions beyond Navier-Stokes order remain open, 
in general, the conceptual basis for hydrodynamics in fluids with elastic collisions seems 
clear and convincing. In the next section, the extent to which a similar analysis applies for 
granular gases is discussed. 



VI. LOCAL HCS AND NAVIER-STOKES EQUATIONS 

The mean free path £ = 1 / ncr^ is a characteristic microscopic length both for fluids with 
elastic and inelastic collisions. Consider a state for which the spatial variations of n, T, and 
U are small on the scale of the mean free path (e.g., iVlnn << 1). Then it is expected 
that the functional dependence of the normal solution on the hydrodynamic fields can be 
made local in space through a Taylor series expansion about the point r and time t for which 
the distribution function is being evaluated. Similarly, since the solution is normal, all time 
derivatives occur only through the fields which obey the balance equations (|5.1|) - (|5.3|) . For 
fluids with elastic collisions, these time derivatives are proportional to gradients of the fluxes 
which in turn are proportional to the small gradients of the fields. As discussed above, this is 
the basis for the dominance of hydrodynamics at long times and long wavelengths. The first 
difference for granular fluids appears at this point in the analysis since the time derivative 
of the temperature is not simply proportional to the gradients but also has a contribution 
from the cooling rate Since this is proportional to (1 — a^) it can be made small 

for sufficiently weak dissipation, although this is not a condition for the Chapman-Enskog 
expansion described below. It is only necessary that ( remains smaller than the decay rates 
for the non-hydrodynamic excitations at large dissipation, a condition for the dominance of 
a closed set of equations including the temperature at late times. 

Let e denote a formal small " uniformity" parameter measuring the small spatial gradients 
in the fields (e.g. a term of order e is of first order in a hydrodynamic gradient, is either 
a product of two first degree hydrodynamic gradients or one second degree hydrodynamic 
gradient). As anticipated above, there is no restriction that the cooling rate ( be of order e; 
only the spatial gradients are being ordered by the uniformity parameter. The distribution 
function, collision operator, and time derivative are given by the representations 

F = F(°) + eF« + ---, J^ = j{o) + eJ« + ..., 9* = ^ + eSf ^ + ■ ■ ■ (6.1) 

The coefficients in the time derivative expansion are identified from the balance equations 
with a similar expansion for Pij, q, and ( generated through their definitions ( [5.7] ) - ( [5. 9] ) as 



functionals of F. The leading term F^^^ is further constrained to have the same moments with 
respect to 1, f^, and v as the full distribution F. To zeroth order in e the macroscopic balance 
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equations become df^n = 0, d^^^u = 0, T ^d^^^T = — C*^°\ and the time derivative in 
the Enskog kinetic equation to this order can be evaluated as 

^(0)^(0) ^ ^^(o)j.Q^pio) ^ l^(o)v^. (v/(°)) (6.2) 

The first equahty follows from the normal form of F^^^ , and the second equality follows from 
dimensional analysis which requires the form 

= n{r, t)vo^^'^(f){V/vo), = 2T(r, t)/m (6.3) 

The dependence on the magnitude of V = v — U{r,t) is due to the requirement that to 
zeroth order in e the distribution function must be isotropic with respect to the peculiar 
velocity. 

The Enskog kinetic equation to zeroth order in e determines the velocity dependence of 
i?(o) 

-C^^^Vy (VF(o)) = JW[FW,fW]. (6.4) 



Comparison of ( |6.3| ) with (O) shows a close similarity with the HCS. In fact, Eq. 



is the same as that for the HCS so the scaling function 0(x) is exactly the same0. Here, 
however, F^'^^ is normal so the density, temperature, and flow velocity of the HCS must 
be replaced by their corresponding exact values for the inhomogeneous state considered, as 
indicated explicitly in (|6.3|) . This is referred to as the local HCS distribution (a misnomer, 
since it is no longer homogeneous). The fluxes can be calculated to this order using only 
the symmetry of this distribution with the results q'-"-* = 0, Pij = p5ij. The pressure p is 
the same as that for the HCS, p = nksT 1 + ^^iina'^g {qu = a) , except now evaluated 

at the true nonuniform temperature and density. The evaluation of the cooling rate 
requires the detailed form for the local HCS, which can be obtained approximately by the 
polynomial expansion method described in Section 3. At low density, i.e. the Boltzmann 
limit, the resulting hydrodynamic equations are the counterpart of the Euler equations for a 
fluid with elastic collisions: first order partial differential equations in space and time. The 
only qualitative difference is the source term —T(^^^ in the temperature equation. At finite 
density, the Enskog equation gives an additional contribution to the cooling rate proportional 
to V ■ U at the next order in the Chapman-Enskog method, which must be retained as well at 
Euler order (hydrodynamic equations to first order in the spatial gradients). The dense fluid 
Euler equations were first obtained from a consistent application of the Chapman-Enskog 
method by Goldstein and Shapiro [|T^ . 

At this point the most significant differences in the application of the Chapman-Enskog 
method to granular gases have been exposed and it is an appropriate point to address some 
misconceptions. Early derivations of hydrodynamics from kinetic theory used the local 



^In Section 3 the Boltzmann equation is considered, while here the Enskog equation is used. 
However, for the HCS and to zeroth order in e the Enskog and Boltzmann collision operators differ 
only by a constant which can be absorbed in (^^^\ 
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Maxwellian for the reference state F^'^\ as is done in the case of elastic colhsions ||46|-^. 
Since the Maxwelhan is only a good approximation for granular fluids at weak dissipation 
the applicability of the resulting hydrodynamics also has this limitation to weak dissipation. 
Similarly, it has been argued that hydrodynamics is a process of relaxation towards equi- 
librium and therefore can apply only if the granular gas is close to the local Maxwellian, 
implying the same limitation. However, the local Maxwellian is not a representation of the 
equilibrium state but rather of the macroscopic non-equilibrium state through its dependence 
on the hydrodynamic fields. Only its velocity dependence is equilibrium-like. Furthermore, 
the analysis here shows that even this velocity dependence is not free to be chosen, neither 
for computational convenience nor to match conceptual bias from experience for fluids with 
elastic collisions. Rather it is determined by the expansion itself. For elastic collisions the 
left side of ( |6.4|) vanishes and the solution is indeed the local Maxwellian as a consequence 
of detailed balance. The physical interpretation is that a wide class of solutions for spatially 
inhomogeneous states evolve in two stages. During a short transient period of the order of 
the mean free time (the kinetic stage), the velocity distribution approaches closely the local 
Maxwelhan and becomes normal as in ( |6.3|) . Subsequently, on a longer time scale the space 
and time dependence of the distribution occurs only through the fields that are governed 
by hydrodynamic equations. A similar interpretation holds for the granular gas. During 
the kinetic stage the velocity distribution rapidly approaches a normal solution close to the 
local HCS, now determined from ( |6.4| ) with the non-zero contribution from the cooling rate. 
The velocity distribution is no longer Maxwellian and the balance equations do not conserve 
energy, but the conceptual basis for this two stage relaxation and the mathematical imple- 
mentation of the Chapman-Enskog method for the normal solution is not compromised by 
these changes. 

It also has been argued that the instabilities of the HCS limit the time scale for validity 
of hydrodynamics. This is a valid question for calculating the linear response of the actual 
HCS to small spatial perturbations. In that case a solution to the kinetic equation is sought 
by expanding about an unstable state. While the derived linear hydrodynamic modes may 
apply through the early stages of instability, eventually the growing perturbation will in- 
validate the description. The Chapman-Enskog method generates a quite different solution, 
which agrees with linear response only when the local HCS is linearized about the HCS. The 
method generates the velocity dependence of a spatially inhomogeneous distribution from 
the Enskog equation through its parameterization by the hydrodynamic fields. Its space 
and time dependence is governed through these fields by the full nonlinear hydrodynamic 
equations. Consequently, it is capable of describing late stage stabilization of linear insta- 
bilities as long as spatial gradients remain relatively small. Examples have been given in 
the previous section for linear instabilities associated with vortex fluctuations that appear 
to evolve nonlinearly to stationary inhomogeneous states. Similar situations occur for fluids 
with elastic collisions, where the hydrodynamic equations describe well a wide range of bi- 
furcation phenomena. In such cases both Euler level nonlinearities and the nonequilibrium 
state dependence of the transport coefficients (temperature and density) can be important, 
and both are incorporated in the Chapman-Enskog method. 

Mathematically, the changes in this method for granular gases arise from the fact that 
the time derivative of the temperature does not vanish to lowest order in e, as it does for a 
gas with elastic collisions. In fact, the reference state incorporates the time dependence 
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of the temperature even for strong dissipation. It has been remarked that the vahdity of the 
Chapman- Enskog expansion is questionable for such a time dependent reference state [Q. 
However, since F^'^^ is normal, it necessarily has the time dependence of all hydrodynamic 
fields even in the case of elastic collisions. The primary difference for granular gases is the 
introduction of a new time scale 1/C''°^ in the reference state. However, there is nothing 
a priori inconsistent with a description of slow spatial decay towards a time dependent 
reference state. Instead, the key question is whether the new time scale is smaller than that 
of other excitations active during the kinetic stage and negligible on the longer time scale. 
If so, the inclusion of the temperature as one of the hydrodynamic fields is justified for these 
long times. Further comment on this is given in Section 8. 

Implementation of the Chapman- Enskog method to the first order in e is now straight- 
forward and has been carried out in detail and without approximation recently for the 
Boltzmann equation and for the Enskog equation [Q; the case of a two component 
mixture is considered in ||^. The constitutive equations for the one component fluid found 
to this order are 



p5ij - r] djU, 



d^U, - -5., V ■ U 



7'5*j V • U, 



(6.5) 



(6.6) 



c 



(1) _ p-p 



■ u. 



(6.7) 



The form of the pressure tensor is the same as that for fluids with elastic collisions, where 
rj (n, T, a) is the shear viscosity and 7 (n, T, a) is the bulk viscosity, depending on the resti- 
tution coefficient as well as the local density and temperature. The heat flux is similar to 
Fourier's law, where k (n, T, a) is the thermal conductivity. However, there is a new trans- 
port coefficient ^{n,T,a) coupling the heat flux to a density gradient, and a coefficient 
{n, T, a) coupling the cooling rate to the divergence of the velocity field. These latter two 
coefficients vanish at a = 1. The transport coefficients are given in terms of solutions to 
inhomogeneous integral equations involving the linearized Enskog operator. Solubility con- 
ditions for the existence of solutions have been proven and approximate solutions in terms of 
polynomial expansions have been obtained for practical purposes, just as for the case of elas- 
tic collisions |^. Consequently the transport coefficients are known as explicit functions of 
the density, temperature, and restitution coefficient. The Chapman-Enskog method places 
no explicit restriction on these parameters, so their domain of applicability with respect to 
density and dissipation is the same as that for the underlying Enskog kinetic equation. In 
particular, there is no restriction to weak dissipation introduced by the method. 

The hydrodynamic equations associated with ( |6.5| ) - ( |6.7| ) are the Navier-Stokes equa- 
tions or Newtonian hydrodynamics (constitutive equations that are linear in the gradients). 
Under some circumstances large gradients occur and more complex constitutive equations 
are required (e.g., shear flow discussed below). The need for more complex constitutive 
equations does not signal a breakdown of hydrodynamics only a failure of the Navier- 
Stokes approximation Although the Chapman-Enskog method can be carried out to 
second order in e (Burnett order), it is likely that failure of the Navier-Stokes approximation 
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signals the need for other methods to construct the normal solution that are not based on a 
small gradient expansion. There are reasons to expect that such non Newtonian transport 
may be more prevalent for granular flows 0. 

The validity of Navier-Stokes hydrodynamics and the dependence of the transport co- 
efficients on the restitution coefficient has been verified in a number of simulations, both 
DSMC and MD, giving good agreement with the predictions from the Chapman-Enskog 
method. The tests at low density based on the Boltzmann equation have been review re- 
cently by Brey and Cubero [Q. The details described there will not be repeated; only a 



summary of the properties studied and the results obtained is given. The Navier-Stokes 
equations described above, with the low density transport coefficients as explicit functions 
of a are considered in reference Using DSMC the evolution of an initial sinusoidal 

transverse component of the velocity field is measured. After a few collisions it is found 
to decay according to shear diffusion, characterized by the shear viscosity. In this way the 
shear viscosity is measured as a function of a for a > 0.6, showing very good agreement 
with the results from the Chapman-Enskog method. Similarly, the longitudinal modes were 
excited by an initial sinusoidal perturbation of the density. The evolution of the density and 
temperature fields, changing an order of magnitude over 30 collision times, is in excellent 
agreement with the theory using transport coefficients evaluated at the chosen a = 0.7 (no 
adjustable parameters). Finally a similar perturbation to excite the diffusive mode of a 
labelled particle allowed measurement of the self-diffusion coefficient in two dimensions; a 



second measurement of the mean square displacement was made for comparison as well [54 



Both results are in excellent agreement with the Chapman-Enskog prediction for a > 0.6 
and with MD results performed for a > 0.7. More recently, the self-diffusion coefficient has 



been compared with MD simulation at finite densities ||55|. The agreement is very good 
for a > 0.6 and na^ < 0.25. Significant deviations are observed at higher densities and 
smaller a. The mean square displacement remains linear in time so that hydrodynamic 
diffusion still occurs, although its description by the Enskog equation is no longer accurate. 
Presumably, this is due to the enhanced shear fluctuations at larger densities and smaller 
a discussed in Section 4. The dependence of the shear viscosity on density and a has been 
compared to results from a DSMC measurement of viscous heating |^. The agreement 
with the Chapman-Enskog prediction is excellent for a > 0.6 and na^ < 0.5 (in this case 
the shear fluctuations are implicitly suppressed by the method). These results show clearly 
the direct veriflcation of hydrodynamics and the quantitative predictions for transport co- 
efficients from kinetic theory for states with small spatial gradients but including strong 
dissipation. 



VII. KINETIC MODELS AND EXACT RESULTS 



In the last Section a precise question was posed: does the Enskog kinetic equation for a 
granular gas imply a corresponding hydrodynamic description on appropriate space and time 
scales. For Newtonian flows and relatively simple nonequilibrium states, the answer appears 
to be "yes". It is possible and of great interest to explore non Newtonian flows as well 
(e.g. steady shear flow), where the Navier-Stokes equations no longer apply. While this is 
possible using DSMC, the Chapman-Enskog method is no longer appropriate to construct the 
normal distribution and associated constitutive equations. This problem has been addressed 
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for the Boltzmann equation with elastic colhsions by considering model kinetic equations 
with structurally simpler collision operators that allow more penetrating analysis. In this 
way a small set of exact normal solutions has been obtained for spatially inhomogeneous 
states far from equilibrium, providing results in semi-quantitative agreement with DSMC 
studies of the underlying Boltzmann equation |^ . Recently, the method of kinetic models 
has been proposed for granular gases based on the Boltzmann and Enskog kinetic equations 
In the former case the Boltzmann equation is rewritten as 

(5t + V- V) /- -(V/) = /(/), (7.1) 



>/'(/) = ^(/)-^^-(V/). (7.2) 



where J{f) is the Boltzmann collision operator defined in (|2.11|) . The velocity derivative 



term on the left side compensates for the time dependence of the temperature in the HCS. 
Consequently, the "new" collision operator J (/)on the right side has properties similar to 
that for elastic collisions: an invariant state and 5 vanishing low velocity moments 

J'(F(°)) =0, J dwi^^ (v) /(/) =0, (l,mv,^ml^2). (7.3) 

Now the usual BGK kinetic model for elastic collisions can be used to represent J'(/) rather 
than J(/), preserving the properties 



J'(/)^-^(/-FW). (7.4) 

where is the local HCS distribution. The collision frequency uj has a space and time 
dependence that occurs only through the density and temperature, and is a free parameter 
of the model. With this choice the model kinetic equation becomes 

(a, + V ■ V) / = -a; (/ - fW) + . (V/) . (7.5) 

A review of some applications of kinetic models for granular gases is given in reference 



Ij]. Such methods are potentially more important for granular gases since the states of 
interest are typically driven by external boundary conditions, posing intractable difficulties 
for analytic studies of the Enskog equation. 

Exact solutions in kinetic theory for spatially inhomogeneous states are exceedingly rare. 
When, furthermore, such a solution is normal and corresponds to a hydrodynamic state far 
from equilibrium an important benchmark is obtained for both conceptual and computa- 
tional issues. In the case of granular gases at low density there are two interesting examples: 



uniform shear flow and planar nonlinear Couette flow [£0[. Both allow controlled dis- 
cussion of nonlinear rheological properties (e.g., shear thinning, normal stresses) that are 
important for a wide class of real granular flows. For uniform shear flow the model results 
predict the components of the pressure tensor as a function of shear rate and a, as well as 
the velocity distribution function itself which is manifestly normal. For nonlinear Couette 
flow the hydrodynamic fields, heat and momentum fiuxes, and velocity distribution function 
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are determined for conditions such that viscous heating dominates coUisional coohng, in- 
cluding large gradients. Explicit expressions for the generalized transport coefficients (e.g., 
viscosity and thermal conductivity) are obtained as nonlinear functions of the shear rate 
and a. Comparison of such non Newtonian transport properties with DSMC results for the 
Boltzmann equation in each case shows surprisingly good agreement, confirming that the 
simplifications of the model do not compromise the physical implications of the Boltzmann 
equation. For the present discussion, these examples demonstrate the existence of a normal 
solution and hydrodynamics under conditions for which the Chapman-Enskog expansion is 
no longer justified. 



VIII. SEPARATION OF SCALES 

It remains to discuss the conditions for dominance of a hydrodynamic description. Simply 
put, the time scale for hydrodynamic processes must be long compared to that required for 
all other ("kinetic") excitations to decay and become negligible. The characteristic longest 
time scale for the kinetic excitations is of the order of the mean free time r^. The analysis is 
simplest for states near the homogeneous state (equilibrium or HCS). Then the amplitudes 
of the hydrodynamic fields do not deviate much from their uniform values, but their spatial 
variation can be complex. For fiuids with elastic collisions the hydrodynamic frequencies, 
LUnik), are functions of the wavevector k (inverse wavelength) such that oJuik) — as 
— *• 0. Consequently, a sufficient condition for the separation of time scales is a restriction 
to sufficiently small k < = (mean free path)~^, or in terms of the frequency wnik) < 
ujq = VQi~^ = collision frequency. Note that this does not require that the mean free 
path should be large compared to the particle size. For example, in dense fiuids the mean 
free path can be smaller than the hard sphere diameter and an accurate (but not Navier- 
Stokes) hydrodynamics applies from macroscopic length scales down to those of the order 
of the particle size. The analysis for inelastic collisions is similar, except that there is the 
additional condition that the cooling rate must be smaller than the collision rate, C < t^O; 
in order to include the temperature as one of the hydrodynamic variables. 

For states near the HCS, the separation of time scales (or frequencies) for a granular gas 
can be studied from the spectrum of the associated linear kinetic equation. There are two 
cases for which precise and unambiguous results have been obtained: the Enskog-Lorentz 
equation for a heavy impurity particle diffusing through a fiuid in its HCS, and a kinetic 
model for the Boltzmann equation. In the first case a mass ratio expansion of the Enskog- 



Lorentz equation leads to a Fokker-Planck equation ||58|j6^ . While there are many qualitative 
differences from the corresponding equation for elastic collisions, there is an exact mapping 
of the case a < to that with a = 1 for all values of the impurity-gas restitution coefficient. 
The spectrum of the Fokker-Planck equation can be obtained exactly in this way. It shows 
a diffusive mode smaller than infinitely many other discrete kinetic modes, bounded away 
from them at a smallest distance of the order of the impurity-gas collision frequency. Fur- 
thermore this distance remains the same independent of the wavelength considered. Thus 
the justification and conceptual basis for hydrodynamic diffusion in this granular "Brow- 
nian motion" is the same as for elastic collisions; it dominates after a few collisions. The 
predictions of this Fokker-Planck description have been confirmed in detail by DSMC for 
the Boltzmann-Lorentz equation and by MD |B3]. 
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Similar conclusions follow for the linearized kinetic model for the Boltzmann equation 
( [r.5|) . As an example, consider an initial state that excites transverse shear excitations (a 
perturbation of the transverse flow field) [Q. The complete spectrum can be determined 
in this well, although it is more complex including branch cuts as well as discrete 

points. There is an isolated shear diffusion mode with dimensionless frequency 



~^ 2 




(8.1) 



where (* = C/ooq is the dimensionless cooling rate. The second line indicates the small 
wavevector limit of {k*,a = 1). This is the unstable shear mode discussed in Section 
4. In this case it is extended to larger wavevectors, beyond the Navier-Stokes form in 
the second line, but there is a restriction on the wavevector to be less than the inverse 
mean free path (more precisely k* = ki < (l + a/tt ), just as for elastic collisions. 
The remainder of the spectrum is at larger frequencies with a smallest separation of the 
order of the collision frequency. Hence shear diffusion dominates after a few collisions. 
The longitudinal excitations can be studied exactly in the same way for this kinetic model, 
with the same conclusion that the hydrodynamic spectrum remains bounded away from the 
kinetic modes and sufficiently small wavevectors. 

In both cases the analysis holds without any restriction to weak dissipation, and the 
cooling rate always remains separated from the kinetic modes. This provides an illustration 
of a contracted description including the temperature as one of the variables, even though 
it is not associated with a conserved quantity. Within the set of hydrodynamic excitations 
there can be a second separation of time scales when k*"^ << (*, i.e. long wavelengths and 
strong dissipation. This is possible since some of the hydrodynamic modes vanish as k*'^ 
while others are coupled to the temperature whose dominant time scale (*~^ is fixed by a. 
Calculations based on the Boltzmann and Enskog equations lead to the necessary condition 
for the temperature to be a hydrodynamic field, (* < Uq. Then hydrodynamics includes 
the two qualitatively different regimes, C* << k*"^ < t^o and k*"^ « (* < uq. Occurrence 
of this separation within the hydrodynamic domain should not be confused with the last 
inequality, which assures the hydrodynamic description. The inclusion of a non-conserved 
degree of freedom, temperature, is not unique to granular gases. Some complex molecular 
gases have internal molecular modes with frequencies comparable to the sound velocity, and 
these modes can be considered as "hydrodynamic" on the dominant long time scale. 

The above analysis refers to states near the HCS, confirming the sufficient conditions of 
hydrodynamic fields whose space and time variation is smooth on the scale of the mean free 
time and mean free path. When the hydrodynamics is non Newtonian, the field variations are 
large and these sufficient conditions are necessarily violated. Nevertheless, it is possible that 
the characteristic hydrodynamic frequencies can still be small compared to the mean free 
path and the separation of time scales can still be valid. This is the case for polymer fiuids, 
where rheological properties imply shear rates large compared to the collision frequency. It 
is difficult to provide a general characterization of non Newtonian hydrodynamics since the 
constitutive equations are not known and may not be universal. However, uniform shear 



fiow is a rare case for which exact results can be obtained |^4[ to illustrate the existence 
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of hydrodynamics for large spatial gradients. An exact, closed set of equations for the 
components of the pressure tensor can be obtained from the Boltzmann equation for the 
special case of Maxwell molecules. The solution shows an exponentially fast approach to 
the hydrodynamic state (uniform shear flow with viscous heating, as obtained from the 
hydrodynamic equations) for times large compared to r (a*) 



Tia =Trr 



2 (l + 9a*2 + 3aV2 + 9a*2)^^^ + 2 (l + 9a*2 - 3a V2 + Qa*^^ 



(8.2) 



where a* = arm is the shear rate times the mean free time, measuring the hydrodynamic 
gradient. It is seen that t (a*) < t^, for all values of a*, so the hydrodynamic description 
applies as expected for t ^ Tm even if a* >> 1. A corresponding approximate analysis for 
hard spheres based on a kinetic model can be extended to inelastic collisions leading 
to similar conclusions. The steady shear flow of a granular gas necessarily entails large a* 
(see comments below) and this has been interpreted as a breakdown of hydrodynamics for 
granular systems However, as the above example illustrates, it is only a failure of 

conditions for the Navier-Stokes form of hydrodynamics, not hydrodynamics more generally 



It is appropriate at this point to remark on the conjecture by Goldhirsch in his recent 
review [Q: "The main thesis of this paper is that granular gases should be considered to be 
mesoscopic in the sense that both the microscopic spatial and temporal scales are typically 
not well separated from the relevant corresponding macroscopic scales and this property of 
granular gases is at the root cause of many and perhaps most of the peculiar properties of 
granular gases." The terminology "not well separated" requires analysis in each case. It is 
not sufficient to make estimates based on dimensional analysis of macroscopic fields in terms 
of characteristic mesoscopic (kinetic) scales. As the example of the last paragraph shows, 
such analysis can lead be misleading (i.e. a* » 1 does not imply the time scales are not 
well separated). The test for the validity of a hydrodynamic description for granular gases 
must be more precise: how do the true hydrodynamic frequencies compare to the true kinetic 
frequencies? For granular gases both the physical states and the forms for these frequencies 
can be more complex than those for normal fiuids. The failure of familiar properties for 
the latter to hold for the former does not necessarily imply a failure of hydrodynamics, and 
there appears to be no precise example of hydrodynamic frequencies exceeding the collision 
frequency. Still, Goldhirsh is correct to warn that granular systems frequently support 
states outside the Navier-Stokes domain and therefore properties based on the latter cannot 
be trusted. 

It is argued in that mean free paths in granular systems can be macroscopic, leading 
to a failure to separate microscopic and macroscopic space scales. The primary example is 
a mean free path defined in the laboratory frame rather than the local rest frame for the 
fiuid element. Such a definition allows arbitrarily large mean free paths for sufficiently large 
convection of the fiuid, and can be anomalously large for both normal and granular fiuids. 
Reasons for rejecting this definition are given in [|66|. In general, when properly formulated 



it appears that the mean free path in a granular medium is qualitatively similar to that in 
a normal fiuid. Of course, in strongly heterogeneous states of dense clusters and large voids 
the average mean free path is no longer a useful concept. 

The title of the review " Granular gases: probing the boundaries of hydrodynamics" 
recognizes an important motivation for studying rapid granular fiow. The hydrodynamics 
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is indeed different from the standard Navier-Stokes form for many reasons, and commonly 
occurring states often sample a different domain of the relevant parameter space than is 
accessible for normal fluids. The single feature of an energy sink allows a wide range of 
phenomena that are not possible otherwise. For example, the steady states of uniform shear 
flow or Couette flow described in Section 7 are possible only because of the balance between 
viscous heating and coUisional cooling. In the steady state, therefore, the temperature, shear 
rate a, and restitution coefficient a are no longer independent variables. In dimensionless 
terms the relationship gives a* oc (1 — o?^ and it is possible to show that the fluid is always 
non Newtonian for any a < 1 [Q. A similar conclusion applies for a granular gas con- 
fined between walls with fixed temperature. The resulting steady state has a nonuniform 
temperature profile resulting from competition between the energy source at the wall and 
collisional cooling in the interior. The profile is fixed by the restitution coefficient and wall 
temperatures. Again, it can be shown that a correct description of the profile cannot be 
obtained from the Navier-Stokes equations due to inherent non Newtonian effects induced 
by the steady state. This is probably a characteristic of most steady states obtained by 
driving the system. Such non Newtonian behavior of simple atomic systems is impossible 
to observe experimentally (except in computer simulations) but seems to be the norm for 
granular gases. This makes the search for complex constitutive equations an interesting 
and challenging research area to be pursued. Rather than abandoning hydrodynamics for 
granular gases, it appears profitable to embrace such a description for new horizons and new 
opportunities. 



IX. BOUNDARY VALUE PROBLEMS 

Most problems of practical interest involve finite geometries and associated boundary 
conditions. Such problems can be addressed directly at the level of the Navier-Stokes hy- 
drodynamics but, as with normal fiuids, become less straightforward for non Newtonian 
fiows. In any case, boundary value problems do not appear to pose any inherent problem 
for kinetic theory or its implementation via DSMC As noted by Grad in his discussion of 
the Boltzmann equation |11|] a normal solution of the form |5.10| can be expected to apply 



only outside certain domains of initial and boundary "slip" (a third slip across a shock layer 
is also discussed). The initial slip layer is the transient period of a few collisions required 
for the kinetic excitations to decay relative to those of hydrodjTiamics, as described above. 
Similarly, there is a domain of the order of a few mean free paths near the boundaries for 
which the solution to the kinetic equation is more complex than that which can be repre- 
sented by the normal form. Away from these slip domains a description in terms of the five 
hydrodynamic fields can be expected. The boundary conditions for these fields necessarily 
require an independent means for extrapolation across the slip domains (modified initial and 
boundary conditions, relative to the actual given conditions). 

For large systems and Navier-Stokes order gradients, the slip conditions are often neg- 
ligible and the interesting physical phenomena resulting from boundary sources are not 



dependent on the mesoscopic details of the boundary layer. Recently, Brey et al. have 
applied the Navier-Stokes hydrodynamic equations derived from the Boltzmann equation as 
described in Section 6 to describe a vibrated granular gas in an external gravitational field. 
It is worth repeating that the transport coefficients and cooling rate are specified functions 
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of the restitution coefficient so there are no adjustable parameters. The system was driven 
at ^ = and unbounded for z > 0. In appropriate dimensionless variables, non-trivial 
temperature and density profiles are obtained independent of the details of the boundary 
conditions (frequency and amplitude of vibration) although a transition to a collisionless 
Knudsen gas eventually occurs at asymptotically large z. These predictions from Navier- 
Stokes hydrodynamics are in very good agreement with corresponding results from both MD 
and DSMC for the Boltzmann equation. Similar good agreement between Navier-Stokes hy- 
drodynamics, MD, DSMC, and experiment has been observed for gravitational flow of grains 
past a fixed wedge [^]. In this case the flow is supersonic leading to formation of a shock 
profile at the tip of the wedge and an expansion fan below the wedge, all well-described by 
the hydrodynamic boundary value problem. 

It is interesting to note that the Navier-Stokes equations are capable of describing some 
of the more exotic behavior exhibited by granular gases, even though structurally they are 
quite similar to the equations for normal fluids. Symmetry breaking is one such phenomenon. 
Consider the above shaken system now partitioned vertically into two domains initially with 
the same density proflles. If a hole in the partition is introduced at a certain height, an 
asymmetry of the density in the two partitions develops depending on the amplitude of the 
vibration and degree of dissipation. The description in terms of a simple effusion process is 



in reasonable agreement with results from experiment and MD |39|. A more detailed Navier 



Stokes analysis of a different symmetry breaking mechanism also is in good agreement with 
MD |7^. In this case a closed container of particles is vibrated at 2; = (no gravity) 



and an exact solution to the Navier-Stokes equations with constant pressure is constructed. 
Next, the system is partitioned normal to z starting at a height z = h. At sufficiently 
low density the hydrodynamic flelds are symmetric with respect to the two sides of the 
partition, but above a critical density an asymmetry in the density occurs (a bifurcation 
in the hydrodynamic solution). The asymmetry (high and low density sides) continues to 
increase with larger density. The hydrodynamic analysis again is in excellent quantitative 
agreement with MD and DSMC results. 

If some dimensions of the system are small, the details of the boundary conditions can 
become more important. For granular systems the boundaries may be rough, e.g. consisting 
of frozen layers of other grains, and the geometry of this roughness (local curvature) also can 
play a role. Coulomb friction is another mechanism for momentum transfer at the boundary. 
The problem of constructing realistic boundary conditions for granular gas hydrodynamics 
has been reviewed and illustrated recently by Jenkins |7T|. As in the case of normal fluids, 
the objective is not to describe in detail the boundary layer near the surface, but rather 
to construct representative conditions for the hydrodynamic equations on the other side of 
the slip. Practical boundary conditions have been found in many cases for the accurate 
application of hydrodynamics for experimental analysis. This is an active and evolving fleld. 



X. DISCUSSION 



Two questions have been addressed here: 1) can kinetic theory provide a valid mesoscopic 
description of rapid flow (fluidized) granular media?, and 2) can a hydrodynamic description 
be formulated and justifled for a macroscopic description? The evidence presented here has 
emphasized idealized conditions and states for which controlled numerical simulations and 
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theoretical approximations can be compared. At this level, it appears that many of the 
questions, concerns, and objections raised regarding the kinetic theory have been removed 
and the domain of validity has been clarified in the past few years. Certainly, kinetic theory 
(Enskog or more sophisticated) appears to be a powerful tool for analysis and predictions 
of rapid flow gas dynamics. Combined with the numerical solution via DSMC, almost any 
boundary or initial value problem for moderate densities can be explored in detail. The 
most practical version of kinetic theory, the Enskog equation, appears to have a somewhat 
more restricted domain in the parameter space of density and restitution coefficient than for 
fluids with elastic collisions, but remains a remarkably rich and accurate basis for analysis 
of moderately dense granular gases. Any observable differences from Enskog provide a new 
opportunity to study correlated collision phenomena whose effects for normal fluids are 
typically small except for the dense fluid. 

The hydrodynamic description for the same class of idealized conditions and states, also 
appears justifled even though the phenomena can be considerably more complex than for 
normal fluids. Nevertheless, this complexity (e.g., instabilities, clustering, rheology) ap- 
pears to be captured quantitatively by a properly formulated hydrodynamics. The required 
separation of time scales, including the temperature fleld, appears justifled even at strong dis- 
sipation based on the few cases for which this question can be addressed quantitatively. For 
weakly inhomogeneous states the Chapman-Enskog method, properly implemented, gives 
both the Navier-Stokes hydrodynamics for a granular gas and the detailed forms for the 
transport coefficients as functions of density and a. Quantitative confirmation of these re- 
sults, as indicated above in many studies, should remove most reservations about the basis 
of hydrodynamics for such states. Non Newtonian hydrodynamics for states in which the 
gradients of hydrodynamic fields are not small is less well-understood, just as for normal 
gases. However, in the latter case such states are rare or unphysical for simple atomic 
systems. In contrast, they appear frequently for granular gases in steady states where the 
gradients are strongly correlated to the coefficient of restitution. Kinetic models and DSMC 
for shear and Couette flow provide instructive examples of non Newtonian hydrodynamics, 
and support for the expectation that hydrodynamics with appropriate constitutive equations 
is applicable to such complex rheological states as well. 

The optimistic view about hydrodynamics for granular media presented here appears at 
odds with the conclusion by Kadanoff in 0: "Can a granular material be described by hy- 
drodynamic equations, most speciflcally those equations which apply to an ordinary fluid? 
It seems to me the answer is "No!"". This conclusion is based on a number of observed 
phenomena: weakly shaken compact material showing composite effects and patterns, com- 
paction (volume reduction) of dense material by tapping, clustering as described in Section 
IV, inelastic collapse (see footnote in Section 2), and one dimensional examples. It is impor- 
tant to note flrst that Kadanoff has in mind Navier-Stokes-like hydrodynamics (universal, 
local, partial differential equations). Indeed, in many of these cases it is difficult to imagine 
a hydrodynamic description of the type requested. These are not the fully fluidized states 
considered here, and it is clear that the presence of compact structures or large voids entail 
heterogenous materials for which Navier-Stokes hydrodynamics was never intended. The 
case for hydrodynamics improves considerably when that term is extended to include more 
general constitutive equations and when the states considered are fully fluidized in two or 
three dimensions. Polymer and colloidal fluids, liquid crystals, and emulsions all admit a 
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hydrodynamic descriptions although more complex than Navier-Stokes, which are neither 
universal nor necessarily local. A single granular material can manifest properties similar 
to many different states of complex fluids, and it is possible that many different forms of 
constitutive equations will be required for an adequate hydrodynamic description of each. A 
single normal material also can have multiple types of macroscopic continuum descriptions 
for different states (fluid, mctastable, elastic). One conclusion here is that granular media 
exhibit a wide range of interesting phenomena for which a Navicr-Stokcs hydrodynamics is 
an accurate and practical tool. More generally, however, there is agreement with Kadanoff 
that granular fluids are in a class of complex materials "with behaviors which are, at this 
moment, not fully understood." Kinetic theory and hydrodynamics (in the broader sense) 
can be expected to provide much of this understanding. 
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